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Abstract. I discuss the multiphase nature of the intracluster medium 
whose neglect can lead to overestimates of the baryon fraction of clus- 
ters by up to a factor of two. The multiphase form of the cooling flow 
equations are derived and reduced to a simple form for a wide class of 
self-similar density distributions. It is shown that steady-state cooling 
flows are not consistent with all possible emissivity profiles which can 
therefore be used as a test of the theory. In combination, they provide 
strong constraints on the mass distribution within the cooling radius. 



1. Introduction 

The multiphase nature of the intracluster medium (icm) in cooling flows was 
demonstrated a decade ago when deprojections of X-ray surface-brightness pro- 
files showed that mass cools and is deposited from the flow in a distributed 
manner, M % r (e.g. Thomas, Fabian &: Nulsen 1987). However, the complex- 
ity of the theory and lack of data with high spatial resolution means that the 
single-phase approximation is still widely adopted. In this paper I attempt to 
present the theory in a palatable form and give examples of its application. 



2. The multiphase cooling flow equations 

The theory of multiphase cooling flows was set out by Nulsen (1986). I rederive 
the equations in a slightly different form here in the hope that they may prove 
more accessible. 

2.1. Derivation of the equations 

We assume an emulsion of density phases which comove with the flow. The 
distribution is described by the volume fraction, f{p,r,t), such that f dp is the 
fractional volume occupied by phases with densities in the range p to p + dp. 
Then / /dp = 1, and the mean density is p = / fpdp. 
Mass conservation gives 

^(p/)+V.(up/) + ^(pp/) = 0, (1) 

where u is the rate of change of position and p is the rate of change of density 
following the flow. The final term in Equation |^ is the equivalent in density 
space of the divergence in velocity space. 
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Integrating over all densities we obtain 

^ + V.u + /3 = 0, (2) 
P 

where 

/3 = 4 lim {ppf). (3) 

This is equivalent to the usual single-phase equation (e.g. Thomas 1988a) except 
that the mass deposition rate is specified in terms of / rather than being a free 
parameter. 

To find how the volume fraction changes with time we use the energy equa- 
tion, 

b IP 'y-l'n?K 

- = —5 + - ?r' 4 

where P is the pressure, re^A is the radiated power per unit volume and 7 = 5/3 
for an ionised plasma. 

Over a wide temperature range appropriate to clusters the cooling function 
can be approximated by a power-law, A oc T", where a ~ 0.5. Then Equation ^ 
can be simplified by moving to a new density variable. Writing 



P 



po(r,t)^x;-V(2-), (5) 



we obtain 



w = {2-a)\ -\w-{2-a) . (6 

\Po iP) 7 P 

If the adiabatic compression term is removed by setting P (x Pq, then the energy 
equation takes a particularly simple form, ib =constant. However, a more useful 
choice is to take po oc p. 

From the final term of Equation |l], we see that at high density when cooling 
is dominant, then ppf ~ constant. This motivates the substitution 

f=^^^w^'--y^^--^g{w,r,t). (7) 
Po 

Then, using Equations ||, ^ and ^ we obtain the following equation for the 
covariant derivate of g (i.e. following the fluid flow): 

9 , fo \Po - a P , , 

- + {3 -a) — + V.u = 0. (8) 

g Po 1 P 

2.2. The form of the density distribution 

In general is a complicated function of position and time. However, we can look 
for solutions in which g has a constant functional form, g = g{w). Only the flrst 
term in Equation ^ depends upon w. Hence we require that g = wdg/dw oc g. 
There are two kinds of solution: 



1. oc exp{—w). This is the most extended distribution which is convec- 
tively stable (it gives P (x Pq). It includes phases of arbitrarily low density. 

2. (7fc oc (1 — w)^~^, < w < 1; k > I. These solutions possess a minimum 
density, p > po- A: = 1 is the least extended, consisting solely of the 
power-law cooling tail. As A; i— >■ oo the solutions resemble goo- 

For other forms of g we must resort to numerical integration to follow their 
evolution. Thomas (1988b) looked at the steady-state evolution of a range of 
distributions with a sharp cut-off at low densities and reached the following 

conclusions: 

• All distributions develop a high-density tail, / they cool. 

• Sufficiently narrow distributions resemble the pure power-law gi by the 
time they begin to be deposited. 

gi and goo bound all reasonable solutions of the cooling flow equations, be 
they self-similar in form or not. I would also argue that plausible formation 
histories for the icm make low values of k more likely than high values. 



3. The Baryon Catastrophe 

Modelling of the icm suggests that the baryon fraction in clusters is of or- 
der O.OQh^^-^ or more (e.g. White ct al. 1993, White & Fabian 1995), where 
h = ifo/100kms~^Mpc~^ is the dimensionless Hubble parameter. However, 
primordial nucleosynthesis limits baryons to a small fraction of the critical den- 
sity, Qb < 0.015 (Walker et al. 1991). If clusters contain a representative 
mixture of material, then it would seem that the Universe must be open with 
density parameter J7 < 0.25 /t~°-^. This has been termed the Baryon Catastro- 
phe. 

One flaw in the above argument is that the modelling is based on single- 
phase models of the icm. We know that there must be a wide range of densities 
at the edge of the cooling flow and it is plausible that the same is true throughout 
the cluster. Moving to a multiphase model can mitigate (but not eliminate) the 
Baryon Catastrophe. 

First let me outline the determination of the gas and total masses of a 
cluster in the single-phase case. The equation of hydrostatic support relates the 
mass of a spherical cluster to the properties of its icm: 

GM{< r) ^ _ 1 dP^ 
r p dr ' 

where G is the gravitational constant and M is the total mass within radius r. 
To determine M we need to know both the density and pressure as a function 
of radius, whereas we usually have only one constraint: the surface brightness 
profile can be deprojected to give the emissivity, ^(r) = n^A. To proceed we 
assume a polytropic relation, P (x ff where 1 < F < | and normalise the 
average, emission- weighted temperature, T^, to some observed value. 
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When moving to a multiphase model the same procedure holds, except that 
we have to average over all density phases. Hence p p = Po f 9 du^, 



and so on. The details of these calculations can be found in Gunn & Thomas 
(1996). Suffice it to say that the form of the solutions is unchanged but the 
gas mass is lowered and the total mass increased over the single-phase case. 
For a = 0.5 the baryon fractions are lower by factors of 0.74 and 0.60 for the 
distributions gi and Qoo, respectively (these factors decrease slightly at lower 
temperatures for which a is lower). 

Finally, let me note that most models ignore complicating factors such as 
magnetic fields and turbulence. Both of these can act to support the gas, thus 
raising the pressure in Equation ^ and lowering the gas fraction even further. It 
is too early to conclude that only low values of J7 are compatible with the data. 

4. Reconstruction of cluster mass profiles 

The self-similar density distributions derived in Section 2 lead to particularly 
simple forms of the steady-state cooling flow equations. I derive these below 
and then show how they can be combined with the emissivity profile to provide 
strong constraints of the mass distribution within the cooling radius. 

4.1. Theory 

Substituting the functional form of gt into Equations |^ and |8| we see that the 
self-similar forms of the cooling flow equations are: 



^ n' 



2A = nlA{To) / w 



(10) 



IP (3 



= 



(11) 



7P {2-a)k 



and 




(12) 



where 




(13) 



To these may be added the equation of hydrostatic support 



VP 



= 



(14) 



P 



where ^> is the gravitational potential. I assume here that the inflow is subsonic — 
this turns out to be a good approximation in all multiphase cooling flow models. 
In a steady-state and spherical symmetry the above equations reduce to 



-T\ — = ^-7> > (16) 

d m r 7 2 — a 



and 

dln-u . 2„ / 1 



dlnr 7 \2 — a 

where 



2 + -S+h^ + A; r, (17) 



^-^■J^ ^^^^ 
is the ratio of the virial to the thermal temperatures, and 

^ ^7 P u k dlnr ^ 

is the ratio of the inflow time to the constant-pressure cooling time. 

Although there appear to be three equations here, the dimensionless ratios 
S and r are the only important variables. The third equation merely acts as a 
scaling (for fixed r, p^~" oc p^~"u). Hence the physics can be captured in just 
two equations: 

dlnE _ ^7 — 1^ T 

d m r 7 2 — a 

and 

= 3 - -[{3 - „) - .(1 - a)lE - (_ + . (21) 

where x = dlnM/dlnr. The g^o equations can be recovered by letting k oc 
and using kr in place of r as the second dimensionless variable. 

The usual way of proceeding is to pick a functional form for the mass profile, 
X, and then to solve for S and r. From this one can generate and emissivity 
profile, C(r), for comparison with the data. Alternatively, we can specify ^(r) 
and (ieier777me the mass profile. Suppose that /3fit = —(l/6)dln^/dlnr is known. 
Then 



dlnr 2 — a + q;7 \{2 — a){2 — a + aj) 

Furthermore this is an eigenvalue problem: requiring that the solution extend 
to r = fixes the outer boundary condition. Hence we can solve for r, S and x- 
We can get a good idea of the behaviour of the solutions by looking at 
the case of constant /Jfif Imposing the physical constraints S > (i.e. a non- 
negative temperature) and x ^ (i-e- mass constant or increasing with radius) 
restricts /3fit to lie in the range 

^ 80 -t- 21 A; 
< /3fit < TT^TTT^ (23) 



2(5 + 3 A;) - ^ ' - 120 + 36 A; 

(in this expression and henceforth I set 7 = 5/3 and a = 0.5 rather than includ- 
ing them explicitly). Thus steep emissivity profiles, Pfn ^ 0.65, are incompatible 
will all steady-state cooling flow models (larger values can occur, however, out- 
side the cooling radius). In addition, The inner value of can be used to 
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constrain k or, if we assume that the virial temperature drops (i.e. S i-^ 0) 
within the cluster core, to measure it. 

Note that the solution for S and r does not depend upon the normalisation 
of the gas and total gravitational masses. If desired, these can be determined by 
fixing the overall temperature and luminosity. The analysis of Section ^ shows 
that the gas density will be slightly lower and the total mass density slightly 
higher than in the equivalent single-phase analysis. 

4.2. Application to A85 

I will now give an example of the application of the theory to the cooling flow 
in the cluster A85. The emissivity profile obtained with the ROSAT HRI (after 
correction for absorption and the instrumental energy response) were kindly 
supplied to me by Clovis Peres. It is shown in 12 arcsec bins in Figure || together 
with a simple broken power-law fit, 

-1 

(24) 

The asymptotic slope of ^ as A; i-^' is very close to the minimum permitted 
for k = 1 which suggests that the solution will require an inner core in the 
mass distribution. This is illustrated in Figure ^. Note that S drops very close 
to zero at r = lOkpc (it tends to a small constant value within this radius). 
The gravitational density profile (i.e. that of the total mass, not just the gas) is 
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Figure 2. The k = 1 cooling flow solution for A85; (a) S (dashed 
line) and r (solid line), (b) the density of the gravitating matter. 

well-fit at radii greater than 20 kpc by a King model, 



Within this radius the density is poorly constrained. Although it appears to rise 
abruptly, only a small change in the slope of ^ would cause it to level off or even 
fall — all we know for sure is that the virial temperature becomes very small. 
Note also that there is only one bin within 20 kpc and this one is most likely to 
be affected by smoothing by the point-spread function, uncertain correction for 
excess absorption, etc. 

The temperature of the gas is approximately constant outside the core ra- 
dius, but drops by a factor of five in to 10 kpc. A temperature decline in the 
centre of clusters is typical of cooling flows observed by ASCA. 

Note that the slope of the mass-deposition profile, r, is close to unity within 
the cooling radius, TcooI ~ 150 kpc. This radius is not a special one for our 
solutions as we have assumed the cooling flow solution holds everywhere. For 
this reason the asymptotic slope of the gravitational density profile at large radii 
should be taken with a pinch of salt. 

The corresponding solution for = cxd is shown in Figure ^. There in no 
core in the gravitational mass profile in this case, with pdark rising as all 
the way into 10 kpc. This is reflected in the temperature profile, however, which 
also rises by a factor of 5 between 150 and 10 kpc. 

I conclude that A85 has a maximum core radius of 120 kpc and that self- 
similar solutions with low values of k provide more plausible temperature gradi- 
ents than those with high values. If the mass profile has a constant-density core 
then A; = 1, as predicted above. 

5. Conclusions 

• Cooling flows are known to be multiphase and the icm of clusters may well 
be multiphase throughout. This can affect mass measurements. 




(25) 
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Figure 3. The k = oo cooling flow solution for A85; (a) S (dashed 
line) and r (solid line), (b) the density of the gravitating matter. 



• I have rederived the cooling flow equations for self-similar density distri- 
butions. Two of these in particular are expected to bound the behaviour 
of all possible flows. 

• The steady-state cooling flow equations are not compatible will all con- 
ceivable emissivity profiles. Thus they can be used as a test of the theory. 

• The solutions provide bounds on M(r) within the cooling radius and given 
k can be used to measure M(r). 

It would be very useful to apply the cooling flow models discussed here to 
galaxies and clusters with well-resolved cooling flows and especially those with 
measured mass-profiles or temperature gradients. 
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